% 


^k 


c 


NBS  TECHNICAL  NOTE  1021 


U.S.  DEPARTMENT  OF  COMMERCE/    National  Bureau  of  Standards 


NATIONAL  BUREAU  OF  STANDARDS 


The  National  Bureau  of  Standards'  was  established  by  an  act  of  Congress  on  March  3,  1901 . 
The  Bureau's  overall  goal  is  to  strengthen  and  advance  the  Nation's  science  and  technology 
and  facilitate  their  effective  application  for  public  benefit.  To  this  end,  the  Bureau  conducts 
research  and  provides:  (1)  a  basis  for  the  Nation's  physical  measurement  system,  (2)  scientific 
and  technological  services  for  industry  and  government,  (3)  a  technical  basis  for  equity  in 
trade,  and  (4)  technical  services  to  promote  public  safety.  The  Bureau's  technical  work  is  per- 
formed by  the  National  Measurement  Laboratory,  the  National  Engineering  Laboratory,  and 
the  Institute  for  Computer  Sciences  and  Technology. 

THE  NATIONAL  MEASUREMENT  LABORATORY  provides  the  national  system  of 
physical  and  chemical  and  materials  measurement;  coordinates  the  system  with  measurement 
systems  of  other  nations  and  furnishes  essential  services  leading  to  accurate  and  uniform 
physical  and  chemical  measurement  throughout  the  Nation's  scientific  community,  industry, 
and  commerce;  conducts  materials  research  leading  to  improved  methods  of  measurement, 
standards,  and  data  on  the  properties  of  materials  needed  by  industry,  commerce,  educational 
institutions,  and  Government;  provides  advisory  and  research  services  to  other  Government 
agencies;  develops,  produces,  and  distributes  Standard  Reference  Materials;  and  provides 
calibration  services.  The  Laboratory  consists  of  the  following  centers: 

Absolute  Physical  Quantities2  —  Radiation  Research  —  Thermodynamics  and 
Molecular  Science  —  Analytical  Chemistry  —  Materials  Science. 

THE  NATIONAL  ENGINEERING  LABORATORY  provides  technology  and  technical  ser- 
vices to  the  public  and  private  sectors  to  address  national  needs  and  to  solve  national 
problems;  conducts  research  in  engineering  and  applied  science  in  support  of  these  efforts; 
builds  and  maintains  competence  in  the  necessary  disciplines  required  to  carry  out  this 
research  and  technical  service;  develops  engineering  data  and  measurement  capabilities; 
provides  engineering  measurement  traceability  services;  develops  test  methods  and  proposes 
engineering  standards  and  code  changes;  develops  and  proposes  new  engineering  practices; 
and  develops  and  improves  mechanisms  to  transfer  results  of  its  research  to  the  ultimate  user. 
The  Laboratory  consists  of  the  following  centers: 

Applied  Mathematics  —  Electronics  and  Electrical  Engineering2  —  Mechanical 
Engineering  and  Process  Technology2  —  Building  Technology  —  Fire  Research  — 
Consumer  Product  Technology  —  Field  Methods. 

THE  INSTITUTE  FOR  COMPUTER  SCIENCES  AND  TECHNOLOGY  conducts 
research  and  provides  scientific  and  technical  services  to  aid  Federal  agencies  in  the  selection, 
acquisition,  application,  and  use  of  computer  technology  to  improve  effectiveness  and 
economy  in  Government  operations  in  accordance  with  Public  Law  89-306  (40  U.S.C.  759), 
relevant  Executive  Orders,  and  other  directives;  carries  out  this  mission  by  managing  the 
Federal  Information  Processing  Standards  Program,  developing  Federal  ADP  standards 
guidelines,  and  managing  Federal  participation  in  ADP  voluntary  standardization  activities; 
provides  scientific  and  technological  advisory  services  and  assistance  to  Federal  agencies;  and 
provides  the  technical  foundation  for  computer-related  policies  of  the  Federal  Government. 
The  Institute  consists  of  the  following  centers: 

Programming  Science  and  Technology  —  Computer  Systems  Engineering. 

'Headquarters  and  Laboratories  at  Gaithersburg,  MD,  unless  otherwise  noted; 

mailing  address  Washington,  DC    20234. 

2Some  divisions  within  the  center  are  located  at  Boulder,  CO    80303. 


Fourier  Transformation  *w  80c,ec 

of  the  Nonlinear  VOR  Model    °°T53 
to  Approximate  Linear  Form      Kf ' 


Dominic  F.  Vecchia 


Statistical  Engineering  Division 
National  Engineering  Laboratory 
National  Bureau  of  Standards 
Boulder,  Colorado  80303 


Prepared  for: 
Department  of  Defense 
Calibration  Coordination  Group 


,<:>.. 


* 


U.S.  DEPARTMENT  OF  COMMERCE,  Philip  M.  Klutznick,  Secretary 

Luther  H.  Hodges,  Jr.,  Deputy  Secretary 

Jordan  J.  Baruch,  Assistant  Secretary  for  Productivity,  Technology  and  Innovation 

NATIONAL  BUREAU  OF  STANDARDS,  Ernest  Ambler,  Director 

Issued  June  1980 


NATIONAL  BUREAU  OF  STANDARDS  TECHNICAL  NOTE  1021 
Nat.  Bur.  Stand.  (U.S.),     Tech.  Note  1021,  32  Pages  (June  1980) 
CODEN:  NBTNAE 


U.S.  GOVERNMENT  PRINTING  OFFICE 
WASHINGTON:  1980 


For  sale  by  the  Superintendent  of  Documents,   U.S.  Government  Printing  Office,     Washington,  D.C.  20402 
Price  $2.00       (Add  25  percent  additional  for  other  than  U.S.  mailing) 


TABLE  OF  CONTENTS 

Page 

1 .  INTRODUCTION 1 

2.  NOTATION 2 

3.  CONTINUOUS   TIME   NONLINEAR  MODEL 3 

3.1  Fourier  Representation   of  General  Signal 4 

3.2  Fourier  Representation   of   VOR  Signal 5 

4.  DISCRETE  TIME    NONLINEAR  MODEL 8 

4.1  Fourier  Coefficients   of   VOR  Signal 8 

4.2  Spectrum   of  VOR  Noise  Model 10 

4.3  Phase  Spectrum  Transformation 11 

5.  LINEAR  MODEL   FOR  PHASE    SPECTRUM 13 

5.1  VOR  Bearing  Angle  Estimator 14 

5.2  Discussion 16 

6.  ROBUSTNESS   OF   ASSUMPTIONS 17 

6. 1  Signal  Offset 17 

6. 2  Signal  Amplitude 18 

7.  ACKNOWLEDGMENT 18 

REFERENCES 19 

APPENDIX 20 


LIST  OF  TABLES 

1.  Bessel  Functions 9 

2.  Squared  Bessel  Functions  and  Weights 15 


1X1 


FOURIER  TRANSFORMATION   OF  THE   NONLINEAR  VOR 
MODEL  TO   APPROXIMATE   LINEAR  FORM 

by 

Dominic  F.  Vecchia 


This  technical  note  describes  a  method  for  transforming  a 
particular  nonlinear  regression  model  to  a  form  which  is  approxi- 
mately linear  in  the  unknown  parameters.   The  technique  involves 
computation  of  the  Fourier  coefficients  for  a  set  of  sample  data 
and  uses  phase  variables  to  estimate  the  parameters.   The  phase 
spectrum  transformation  is  employed  to  obtain  bearing  angle  esti- 
mates for  a  model  associated  with  the  Very-High-Frequency  Omni- 
Directional  Range  (VOR)  aircraft  navigation  system.   The  trans- 
formation provides  a  model  linear  in  relevant  phase  parameters. 
Thus,  estimation  of  VOR  bearing  angle  utilizes  existing 
statistical  theory.   Finally,  it  is  shown  that  certain  generaliza- 
tions of  the  VOR  model  also  are  reduced  to  approximate  linear  form 
by  the  phase  spectrum  transformation. 

Key  Words:  Fourier  coefficients;  linear  model;  nonlinear  model; 
phase  spectrum  transformation;  spectrum;  VOR  aircraft  navigation 
system;  white  noise. 


1.   INTRODUCTION 


There  are  many  different  reasons  for  making  a  transformation  of  variables 
in  the  statistical  analysis  of  data.   This  technical  note  discusses  an  unusual 
type  of  transformation  useful  in  connection  with  a  particular  nonlinear 
regression  model  for  audiofrequency  signals  from  the  Very-High-Frequency 
Omni-Directional  Range  (VOR)  air  navigation  system.   The  model,  which  is 
considered  in  a  more  general  form  than  the  VOR  requires,  is  intrinsically 
nonlinear  in  the  unknown  parameters.   By  intrinsically  nonlinear  we  mean  that 
a  single  observation  cannot  be  transformed  into  linear  form.   For  example, 
consider  the  two  models 

Y  =  exp(B]_  +  62x  +  £) 

Y  =  [e1/(S1-62>  ]  [exp(-B2x)  -  exp(-BlX)]  +  e  , 

where  3}  and  $2  are  unknown  parameters,  x  is  an  independent  variable,  and  e 
is  a  random  error  term.   Both  models  are  nonlinear  in  &\   and  32,  but  the  first 
is  intrinsically  linear  because  the  transformed  variable  lnY  is  linear  in  $1 
and  32«   However,  the  second  model  is  intrinsically  nonlinear  because  it  is 
impossible  to  convert  the  model  into  a  form  linear  in  the  parameters.   For  a 
discussion  of  these  concepts  see  reference  [1]. 

Usually,  it  is  not  useful  to  transform  a  model  of  the  second  type  because 
it  remains  nonlinear  whatever  transformation  is  applied.   The  transformation 
introduced  in  this  paper  is  unusual  because  it  involves  computation  of  the 


Fourier  coefficients  of  the  sample  data  and  uses  phase  variables  to  estimate 
parameters.   For  this  reason  the  procedure  to  be  described  is  called  the  phase 
spectrum  transformation.   The  method  will  be  demonstrated  for  the  model 
specific  to  the  VOR  air  navigation  system. 

The  VOR  is  a  fundamental  component  of  the  present-day  air  navigation 
system.   A  feature  of  the  VOR  system  which  provides  much  versatility  for 
defining  controlled  airways  is  that  the  facility  emits  an  infinite  number  of 
radial  courses  providing  aircraft  bearing  information.   This  information  is 
contained  in  the  phase  angles  of  two  30  Hz  audiofrequency  signals.   The  first 
has  a  constant  phase  at  all  points  around  a  VOR  station  and  is  called  the 
reference  signal.   The  other,  called  the  variable  signal,  has  a  phase  equal  to 
the  bearing  angle  to  (or  from)  the  VOR  transmitter.   In  the  aircraft,  bearing 
information  is  determined  by  measuring  the  phase  difference  between  the  two 
component  signals. 

The  importance  of  the  accuracy  of  bearing  angle  estimation  devices  can- 
not, of  course,  be  overstated.   At  present,  measurement  accuracy  for  VOR  test 
instruments  depends  upon  calibration  with  commercial  equipment  designed  for 
that  purpose.   As  system  requirements  become  more  severe  because  of  increasing 
traffic  in  the  air  lanes,  it  is  clear  that  both  the  accuracy  and  precision  of 
present  VOR  calibration  equipment  will  require  additional  scrutiny.   Hope- 
fully, this  will  increase  the  safety  and  efficiency  of  aircraft  operations. 
For  a  general  discussion  of  the  VOR  system,  see  [2]. 

This  paper  presents  a  statistical  technique  for  estimation  of  VOR  bearing 
angle  and  gives  the  corresponding  precision  of  the  estimated  angle.   The 
general  method  is  based  on  regression  analysis  of  samples  taken  by  a  sample- 
and-hold  amplifier  and  an  analog-to-digital  converter.   The  method  provides  a 
linear  model  in  relevant  phase  parameters.   Thus,  the  bearing  angle  estimation 
utilizes  existing  statistical  theory. 

In  section  3  of  this  technical  note,  the  nonlinear  regression  model  is 
represented  in  continuous  time.   Fourier  coefficients  are  obtained  for  the 
noise-free  signal,  and  results  for  the  special  case  of  the  VOR  signal  are 
stated.   The  results  for  the  VOR  application  are  extended  without  proof  to  the 
discrete  time  sample  model  in  section  4.   The  spectrum  for  the  VOR  model  with 
noise  is  derived  and  the  phase  spectrum  transformation  is  defined.   In  section 
5  the  approximate  linear  model  for  transformed  variables  is  used  to  estimate 
unknown  parameters  and  the  usefulness  of  the  transformed  model  is  discussed. 
Section  6  is  a  limited  discussion  concerning  the  properties  of  estimators  if 
some  assumptions  specific  to  the  VOR  model  are  invalid. 


2.   NOTATION 


For  X,  a  random  variable  with  probability  density  f (x) ,  we  denote  the 
mean  and  variance  of  X  by  E[X]  and  Var[X] ,  and  the  covariance  between  X  and  a 
random  variable  Y  is  denoted  by  Cov[X,Y]. 

Vectors  and  matrices  are  denoted  by  underlined  letters,  for  example,  _9 
and  _V.   If  9  denotes  a  vector,  then  9'  will  denote  the  transpose  of  _9.   An 
estimator  of  9  will  be  denoted  by  9. 


3.   CONTINUOUS  TIME  NONLINEAR  MODEL 


The  nonlinear  model  considered  in  this  report  will  be  represented,  ini- 
tially, as  a  continuous  function  of  time.   In  a  later  section  the  mathematical 
results  obtained  for  the  continuous  time  model  are  extended  to  the  case  where 
the  data  are  equally  spaced  observations  from  the  continuous  signal. 

To  represent  the  deterministic  component  of  the  model  requires  two  peri- 
odic functions  described  below  in  (3.0.1).   These  functions  are  added  to 
obtain  the  expected  (ideal)  value  of  the  output  signal  in  a  nonlinear  regres- 
sion model.   The  two  component  functions  are 

v(t;M)  =  «!  cosUTTf^t+fi)  +  <b{] 
and  (3.0.1) 

s(t;S,J>)  =  a2  cos[2TTf2(t+6)  +  <j>2  +  6  sin[2irf ^t+6)  +  <j>3]  ]  . 

The  waveform  generated  by  the  sum  of  v(t;S,<J>)  and  s(t;<5,<f>)>  with  some  param- 
eters assumed  known,  may  be  used  to  represent  the  ideal  audiofrequency  signal 
for  the  VOR  aircraft  navigation  system.   In  this  context,  v(t;6,_^)  is  called 
the  variable  phase  signal,  and  s(t;<5,_<J>)  is  the  frequency  modulated  subcarrier 
signal.   The  frequency  modulating  sinusoid,  contained  in  the  argument  of 
s(t;6,j0,  is  called  the  reference  phase  signal.   Equations  for  the  component 
signals  have  been  presented  in  a  more  general  form  than  the  VOR  application 
requires.   However,  the  above  terminology  is  used  throughout  this  paper.   Fol- 
lowing are  descriptions  of  the  model  parameters: 

o^  =  amplitude  of  variable  phase  signal; 

a2  =  amplitude  of  subcarrier  signal; 

ll   =  variable  (and  reference)  signal  frequency; 

f2  =  subcarrier  frequency; 

3  =  modulation  index; 

<5  =  arbitrary  fixed  time  offset; 

<J>1  =  phase  angle  of  variable  phase  signal; 

<J>2  -   phase  angle  of  subcarrier  signal; 

<j>3  =  phase  angle  of  reference  phase  signal. 

The  fixed  time  offset,  6,  is  included  in  (3.0.1)  because  the  output  signal 
will  be  observed  and  sampled  from  an  unknown  starting  point  in  the  waveform. 
We  cannot,  in  general,  be  assured  that  observation  of  the  signal  begins  at  a 
zero  crossing  on  the  time  axis. 

Realistically,  measurement  of  the  composite  signal  involves  random 
measurement  error  in  some  form.   In  this  paper  the  random  error  process,  e(t), 
is  assumed  to  be  additive  white  noise  [3],  and  the  resulting  process  Y(t)  can 


be  represented  by 

Y(t)  =  V  +   v(t;6,_£)  +  s(t;M)  +  e(t)  ,  (3.0.2) 

where  P  is  a  fixed  but  unknown  offset.   Unless  otherwise  specified,  e(t)  is 
not  assumed  to  be  Gaussian  white  noise.   However,  for  the  distributional 
result  obtained  in  the  appendix,  we  require  normality  and  independence  of  the 
discrete  time  error  series  to  be  as  described  in  subsection  4.1. 

In  subsection  3.2,  a  form  of  the  model  specific  to  VOR  navigation  system 
is  discussed.   For  this  application  some  of  the  parameters  in  the  general  form 
of  the  model  are  assumed  to  be  fixed,  known  values.   On  this  basis  a 
statistical  method  for  VOR  bearing  angle  estimation  will  be  derived.   Because 
the  relevant  angle  for  the  VOR  application  is  C  4*1 — 4>3 )  »  it  is  sufficient  to 
consider  a  reparameterized  form  of  the  model  where  we  define  new  parameters  _8 
by 

©1  =  4>l  —  <t>3 
92  =  2TTfL6  +  <f>3 
83  =  27rf26  +  <j>2  . 
For  this  parameterization  the  general  model  becomes 

Y(t)  =  u  +  v(t;_6)  +  s(t;_6)  +  e(t)  (3.0.3) 

where 

v(t;_6)  =  <xi    cos[2Trf1t  +  9j  +  02] 
and 

s(t;_9)  =  ot2  cos[2irf2t  +  83  +  3  sin[2irf1t  +  82]  ]  . 

In  the  following  section,  we  obtain  the  Fourier  sine  and  cosine  transforms  of 
E[Y(t)]=y+v(t;_8)+s(t;j3)  under  the  assumption  that  f2=mfi«   Hence,  E[Y(t)]  is 
periodic  with  period  (1/fj).   Utilizing  the  general  result,  the  specific 
transform  for  the  VOR  signal  is  determined  in  subsection  3.2. 

3.1   Fourier  Representation  of  General  Model 

Let  y(t;_8)  denote  the  expected  value  of  Y(t).   The  deterministic  function 
y(t;_8)  is  given  by 

y(t;_8)  =  v  +   v(t;_8)  +  s(t;_8). 

Suppose  that  the  frequencies  f^  and  f2  in  the  definition  of  v(t;_8)  and  s(t;_8) 
are  such  that  f2=mf^  for  some  positive  integer  m.   Thus,  y(t;_6)  is  periodic 
with  period  (1/f^).   Under  this  assumption,  the  Fourier  coefficients  of  y(t;_8) 
can  be  obtained  from  the  real  and  imaginary  parts  of  the  complex  integral 


l/2f1 
sk  =  2fj   /   y(t;_8)  exp[i2uf1kt]  dt,  k=0,±l,±2,  ... 
-l/2fj 

Let  ak  and  b^.  denote  the  real  and  imaginary  parts  of  sk,  so  that  sk  =  a^+ib^. 
That  is,  ak  and  b^  are  the  Fourier  cosine  and  sine  transforms,  respectively, 
of  y(t;_8).   In  the  appendix,  it  is  shown  that  the  Fourier  coefficients  ak  and 
bk  for  k=0,l,...,  are 


and 


r    V  ,    k=0 

ak   =  "S      al    cos  [8^+62]    +   a2   ak(s)    ,    k=l 

^    a2   ak(s)  ,    k>2 


'0  ,    k=0 

-o^    sin[  8^+62]    -   a2   b^Cs)    ,    k=l 


(3.1.1) 


bk  =< 


v. 


-a2   t^Cs) 


,    k>2 


where,    for  k>l, 


ak(s)   =  Jk_m(B)cos[83  +  (k-m)82]   +  J_k_m(  B)cos  [  83  -   (k+m)82] 


and 


b^s)    =  Jk_n(B)sin[83  +   (k-m)82]    +  J-k_m( B)sin[ 83   -   (k+m)82] 


The  notation  Jn(z)  denotes  a  Bessel  function  of  the  first  kind  [4], 

Equations  (3.1.1)  represent  the  Fourier  coefficients  for  the  mean  value 
function  of  the  general  continuous  time  model  where  the  subcarrier  frequency 
f2  is  an  integer  multiple  of  fj,  the  frequency  of  the  variable  phase  signal. 
For  the  VOR  model,  where  some  parameters  in  the  general  model  are  assumed 
known,  we  will  see  that  the  Fourier  coefficients  can  be  greatly  simplified. 
The  form  of  the  coefficients  for  this  special  case  will  suggest  a  method  for 
estimating  the  unknown  parameters. 


3.2   Fourier  Representation  of  VOR  Signal 

The  VOR  audio  frequency  waveform  consists  of  a  30  Hz  variable  phase 
signal  linearly  added  to  a  frequency  modulated  9960  Hz  subcarrier  signal.   The 
modulation  index  for  the  reference  phase  signal  is  assumed  to  be  fixed  and 
known,  as  are  the  amplitudes  of  each  signal.   Specifically,  parameter  values 
assumed  to  be  known  are: 


u  =  0 
<*!    =   21/2 
a2   =   2^2 
fx   =   30 

f2   =  9960  (so  ra  =   f2/fi   -  332) 
3     =  16     . 
These    specifications   define   the  VOR  model 

Y(t)    =  y(t;_6)  +  e(t) 
where 

y(t;_9)  =  21/2{cos[2iT30t  +   Bl  +   92] 

+  cos[2Tr9960t  +  G3  +  16sin[2Tr30t  +  92]]}. 

Noting  that  y(t;9)  is  nonlinear  in  the  unknown  parameters  9"=[ 8^, 92, 93] ,  we 
will  show  in  section  4  that  the  spectrum  of  Y(t)  can  be  used  to  transform  a 
set  of  sample  data  to  new  observations  satisfying  a  model  approximately 
linear  in  8.   The  transformation  to  linearity  will  depend  on  the  simplified 
form  of  the  Fourier  coefficients  for  y(t;9)  when  known  values  of  parameters  in 
the  VOR  model  are  substituted  in  the  general  equations  (3.1.1). 

Substituting  known  values  in  the  expressions  for  ak(s)  and  ty^s),  we 
obtain  for  k>l, 

ak(s)  =  Jk_332(16)cos[93+(k-332)e2]  +  J_k_332( 16)cos [ 93-(k+332) 92] 

and 

bk(s)  =  Jk_332(16)sin[93+(k-332)92]  +  J_k_332( 16) cos [ 83-(k+332) 92] . 

We  need  the  following  results  for  Bessel  functions  of  the  first  kind. 

Lemma  3.1;  [4,  page  358]   For  integer  n,  Jn(z)  satisfies 

J_n(z)  =  (-l)n  Jn(z)  . 

Lemma  3.2;  [4,  page  365]   For  fixed  z,  as  n-*-00  through  real  positive  values, 

Jn(z)  S  (2Tin)-1/2(ez/2n)n  . 


From  lemmas  3.1  and  3.2,  it  follows  that 

J_k_332(16)  =  (-l)k+332Jk+332(16) 

=   (-l)k+332[2ir(k+332)]-1/2[16e/2(k+332)]k+332 

=  0 


Clearly,    the   value   of   Jk+332(16)    i-s    immeasurably    small.      Note,    for   example, 
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that  if  k=l,  Jk+332(16)  ~   •02( * 065)    •   Similarly,  for  small  k,  we  have  that 
Jk-332  -  °- 

Using  the  above  results,  we  get 


and 


ak(s)    = 


0  ,   k=l 

Jk_332(16)cos[e3  +  (k-332)92]    ,   k^2 


T>k(s)    = 


0  ,   k=l 

Jk_332(16)sin[63  +  (k-332)92]    ,    k>2 


Substituting   in   the   general   expressions   (3.1.1),    the   approximate  Fourier 
coefficients   for  y(t;9)   in   the  VOR  model   become 


and 


r  o 


ak  M 


1/2 


2i/z  cos[e1  +  e2] 


k=0 


,  k=l 


^_   2i/z  Jk_332(16)cos[63  +   (k-332)G2]     ,    k>2 


r 


1/2 


bk    2<    -2W/   sintSj   +   92] 


,    k=0 
,    k-1 


v 


-21/2  Jk_332(16)sin[63  +  (k-332)G2]    ,    k>2. 


Because   omitted   terms    are   negligible,    in   the   following   sections   we    consider 
the  Fourier    coefficients    to   be   exact. 


DISCRETE  TIME   NONLINEAR  MODEL 


The   general  model  and    corresponding  Fourier   transforms    introduced   in 
section   3  will   facilitate   a   later   discussion  about   errors    in  assumptions    for 
the  VOR  model.      Because   the   approximations   discussed   in   the   previous    section 
depend   on   the   particular   values    of    some   parameters    in   the   general  model,    the 
phase   spectrum   transformation  will   be   developed   only    for   the  VOR 
specifications.       In   subsection   4.1   we    consider   the   discrete    time    analogs    of 
the   VOR  model   and    corresponding  Fourier    coefficients,    since   a   digital  phase 
estimation   technique    is   desired. 


4.1      Fourier  Coefficients   of  VOR  Signal 

Let   Y-;,    j=l,...,N,    be  N  equally    spaced   observations   from  one   period   of 
the    continuous    time   VOR  series.      For    convenience,    N   is   assumed   to   be   even   in 
the   results    to   follow.      The    sample  model   for   VOR  applications    can   be  written 

Yj  =   yjd)   +  ej  I 

_  >     J»k=l,...,N 

E[ej]=0;    Var[e-j]=a2;   E[ejek]=0   if   j*k     J 

where 

yj(_9)    =  21/2{cos[2tt(j-1)/N  +    d1  +   62] 

+   cos[2tt332(j-1)/N  +    93  +   16sin[2ir(  j-l)/N  +    62]  ]}. 

2 
The   ei's   denote  un correlated   random  error   terms   with  unknown   variance    a   . 

The  Fourier   coefficients   (sk  =  a^+ib^)    for  yj(_9)  are  given  by 

N 
sk   =    (2/N)    I  yi(_6)   exp[i2TTk(j-l)/N]    . 

J-l 

The  ak  and  b^  follow  directly  from  the  continuous  time  model.   Excluding  the 
coefficients  for  k=0,  which  are  not  useful  to  estimate  _9,  equations  for  ak  and 
bk  are 


ak 


21/2    cos [6!  +   82]  ,   k=l 

21/2  Jk_332(16)cos[e3  +  (k-332)92],   k=2,3, . . .  ,N/2, 


and  (4.1.1) 

f  -21/2  sinlSj  +  e2]  ,   k=l 

.  -21/2  Jk_332(16)sin[e3  +  (k-332)62],   k=2,3, . . . ,N/2-l    . 


K 


Recall  that  these  coefficients  can  reasonably  be  considered  exact  expressions 
because  omitted  terms  are  negligible.   The  interesting  feature  of  the 
equations  for  ak  and  bj^  is  the  form  of  the  phase  of  the  kth  harmonic.   If  we 
let  q  represent  the  phase  at  a  chosen  harmonic,  then  the  basic  equation  for 
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q  is  tan(q)=-b/a.   For  convenience,  we  have  dropped  the  subscripts  on  q,  b, 
and  a.   Because  Arctan(-b/a)  gives  the  same  value  for  -b  and  -a  as  for  b  and 
a,  the  full  solution  for  q  in  the  interval  (-t,tt]   is  the  following 
[5,  page  12]  : 


r 


i*  = 


< 


Arctan(-b/a)   ,  a>0 
Arctan(-b/a)-ir,  a<0,  b>0 
Arctan(-b/a)+iT,  a<0,  b<0 


(4.1.2) 


-tt/2 

,  a=0,  b>0 

tt/2 

,  a=0,  b<0 

arbitrary 

,  a=0,  b=0 

The  notation  Arctan(x)  is  used  to  denote  the  principal  value,  so  that 
-ir/2<Ar ctan(x)<iT/2.   If  arctan(x)  denotes  any  angle  whose  tangent  is  x,  then 
it  follows  that 


arctan(-b/a)  =  q*  +  j2ir 


(4.1.3) 


where  j  is  an  arbitrary  integer.  From  (4.1.3)  and  the  expressions  for  a^  and 
b^,  it  then  follows  that  there  exist  integers  Y^  such  that  qk(_9)=q^*+Yic2'T  is 
given  by 


qk<_8)    = 


»1   +    92 


,    k=l 


(k-332)92  +   93    ,    |k-332|<K      , 


(4.1.4) 


where  K  is   a    constant    chosen   to  assure   that   Jk-332^-'-^)    ^s   non-negligible.      A 
zero  value   for   the  Bessel  function  leads    to  an  arbitrary   phase   because 
ak=kk=0*      Table   1   lists   the  values   of  Jn(16)    for  n=0,...,24.      In  a   later 
section   it  will   be   shown    that   a  value   of   K^IO   is   sufficient    for   the  proposed 
estimation   of    9. 


Table    1.      Bessel   Functions 


n 


Jn(16) 


n 


Jn(16) 


0 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

12 


-.1748990739 

.0903971756 

.1861987209 

-.0438474954 

-.2026415317 

-.0574732704 

.1667207377 

.1825138237 

-.0070211419 

-.1895349656 

-.2062056944 

-.0682221523 

.1124002349 


13 

.2368225047 

14 

.2724363352 

15 

.2399410820 

16 

.1774531934 

17 

.1149653049 

18 

.0668480795 

19 

.0354428740 

20 

.0173287462 

21 

.0078789915 

22 

.0033536066 

23 

.0013434266 

24 

.0005087450 

The  multiples  of  2ir  indexed  by  Y  in  the  expression  for  q^(_9)  are 
necessary  to  adjust  the  qk*  from  the  interval  (  —  tt,  tt]  to  the  interval  (-00,00). 
Because  the  9's  will  represent  unknown  parameters,  the  Y's  are  not  known  in 
general.   If  we  assume,  however,  that  02^0,  it  is  clear  from  (4.1.4)  that  the 
Y's  must  satisfy 

qj*  +  Yj2tt  <  qj*  +  Y±2tt  ,  j<i  ;  i*l*j   , 

so  (4.1.5) 

(Yj   "    Yi)2TT    <  q±*   -    qj*    ,    j<i    ;    i*l#j      . 

This    implies    that    for    | k— 332 | <K,    the   Y^'s    can   be   uniquely   determined   if,    for 
example,    the    following    conditions   are   assumed: 

-ir<91+92<TT 

62>0 

-TT<e3<TT      . 

With   these    constraints    Yi=Y332=0   and  successive   values    of    Yk   near  k=332 
are   determined   from   (4.1.5).      For   the   extension   of    the   results    to   the  VOR 
signal  with   error  we   will   assume    that   the   Y's   are  known. 

4.2      Spectrum   of  VOR  Noise  Model 

In   the   previous    subsection  Fourier    coefficients   were   obtained   for   the 
deterministic   component    of    the  VOR  process.      For   the   nonlinear  VOR  regression 
model: 

Yj   =   yj(_e_)   +  ej  1 

2  >     J,k=l,...,N 

E[e-j]=0;    Var[ej]=a   ;    E[ejek]=0    if   j*k      J 

the  Fourier    coefficients    of    y  i(®) ,    represented   by    s^a^+ib^,    were   shown   to 
have    phase   values    linear    in  _6.      In   this    section  we    prove    that    random  phase 
variables   derived  from  Fourier   transformation  of  Y-j,    j=l , . . .  ,N  are  appropri- 
ately  represented   by   a   regression  model   linear   in    9. 

Letting   the    random  variables  S^=A^+iB}c,    k=l , . . .  , N/2-1 ,    represent   the 
Fourier   coefficients   of  Y-j,    we   have 

N 
Sk   =    (2/N)    I  Yj   exp[i2irk(j-l)/N] 
j=l 

N 
=    (2/N)    I    (yi(_9)+ei)exp[i2TTk(j-l)/NJ 

j-l 

N 
=    sk  +   (2/N)    I  ej   exp[i2irk(j-l)/N] 

j=l 
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where   sic=aic+ibjc  are   the  Fourier    coefficients   for   y-j(_9)   given   in  (4.1.1).      If 
we   let    the   transform  of    the   random  error   sequence   be   denoted   by 

N 
gk  +  ihk  =    (2/N)    I  eA   exp[i2irk(j-l)/N]      , 
j=l 

it   is   well  known    [6,    page    110]    that    the   random  variables   gk   and  hk, 
k=l , . . . , N/2-1    are  mutually   un correlated  and 


E[gk]    =  E[hk]    =   0 
Var[gk]    =  Var[hk]    =   (2/N)</ 


k=l,... , N/2-1 


It  therefore  follows  that  the  Ak  and  Bk  are  uncorrelated  and  satisfy  regres- 
sion models 


Ak   =  ak  +  §k 
Bk  =   bk  +  hk 


k=l,... , N/2-1 


Linearity   of   ar ctan(-bk/ak)    in  _6   suggests    that   we    consider   the   phase   spectrum 
of  Yj    to   estimate  _9.      In   the   next    section   it    is    shown   that    the   expected  values 
of   phase   variables   are   approximately    linear   in    8. 

4.3      Phase  Spectrum  Transformation 

The   definition   of   phase   random  variables   parallels    the   description   of    the 
phase   qk(_9)   for   the   deterministic   component   y ^ (_9 ) •      Phase   variables   will   be 
denoted   by   Qk   and   initially   are   defined  using  principal  values    in   the   interval 
(  — 7T ,  tt]  .      Define 


0^.*   =  Arctan[-Bk/Ak]    ,    k=l , . . .  ,N/2-l    . 


(4.3.1) 


In  the  appendix  the  distribution  of  the  (\* ,  k=l , . . . , N/2-1  is  determined  when 
the  errors  are  Gaussian.   The  expected  value  of  0^.*  is  not  obtained,  but  the 
complexity  of  the  distribution  illustrates  the  usefulness  of  approximate 
moments  of  Q^*  which  result  from  a  suitable  propagation  of  errors  formula. 

To  conclude  this  subsection  we  obtain  approximate  formulas  for  the  mean 
and  variance  of  Qk*«   These  results  are  the  basis  for  the  linear  models  used 
to  estimate  _9  in  Section  5.   As  defined,  Ak  and  Bk  appearing  in  (4.3.1) 
satisfy: 


E[Ak]  =  ak  ;  E[Bk]  =  bk 
Var[Ak]  =  Var[Bk]  =  (2/N)a' 
Cov[Ak,Bj]  =  0   all  j,k 


~\ 


V   k,j=l,...,N/2-l 
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For   values   of   k   such   that   a^O^b^.,    it    can   be   shown    [7,    page   333]    that, to   order 
N      ,   Qk*  has   approximate  mean  and  variance   given   by: 

E[Qk*]    =  Arctan[-bk/ak] 

Var[Qk*]    3   (2/N)rk_2a2 

where 


(4.3.2) 


2  2  2 

rk     =  ak     +  bk 


2  ,      k=l 

2J2k_332(16),       |k-332|<K  . 


The   value    of   K   is    chosen   to   assure    that   Jk-332(16)    is   nonzero.      Note   that 
E[Qk*]-qk*,    where   the    solution   for   qk*   in   the   interval   (  —  it  ,  tt]      was   given   in 
(4.1.2).      Based   on   the   discussion   following    (4.1.2)   we    can   define   Qk=Qk*+Yk2T 
where   the   Yk's   are   (generally)   unknown   integers    such   that 

E[Qk]    =   qk(_9) 
Var[Qk]    =   (2/N)rk_2a2 
where 


(4.3.3) 


The  method  of  estimating  _9  to  be  outlined  in  Section  5  is  based  on  the 
above  results.   To  the  chosen  degree  of  approximation,  the  important  features 
are: 

1.  E[Qk]  is  linear  in  _6  for  all  permissible  k. 

2.  Var[Qk]  is  proportional  to  a2  with  known  constant  of  proportionality. 

3.  For  all  permissible  k,  the  Qk's  are  mutually  uncorrelated  (assuming 
white  noise  errors  in  the  original  model). 

Because  the  linear  model  of  Section  5  requires  that  Qk  be  an  observable  random 
variable,  it  will  be  necessary  to  assume  that  the  Yk's  are  known  integers. 
The  Yk's  are  used  to  adjust  the  computed  values  of  the  arctan  function  to 
satisfy  inequalities  implied  by  (4.3.3).   Though  it  is  not  obvious  that  this 
correction  can  be  accomplished  with  the  C^,  which  are  subject  to  error, 
computer  simulations  indicate  that  the  adjustment  is  possible  if  the 
measurement  error  variance,  a2,  is  small.   Values  used  for  a2  are  thought  to 
be  representative  of  measurement  precision  for  a  new  system  designed  to  obtain 
sample  values  Y^,  j=l,...,N. 
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5.   LINEAR  MODEL  FOR  PHASE  SPECTRUM 

In  the  previous  section  approximate  formulas  for  the  means  of  phase 
random  variables  for  the  nonlinear  VOR  model  were  shown  to  be  linear  in  the 
unknown  parameters  _9.   Corresponding  variance  approximations  are  unequal  at 
the  harmonics  and  proportional  to  a  ,  but  do  not  depend  on  other  unknown 
parameters.   The  additional  observation  that  phase  random  variables  are 
uncorrelated  suggests  that  9  and  a     may  be  estimated  using  a  linear  model  in 
with  known  error  covariances.   For  the  description  to  follow  the  reader  is 
reminded  that  expressions  for  the  mean  and  variance  of  Q^  are  not  stated  as 
approximations. 

Consider  the  n=2+2K  equations 

r     Q1  +   82  +  e1  ,  k=l 

Qk  =  \ 

I    (k-332)92  +  63  +  ek   ,  k=332-K, ... ,331,332,333, ... ,332+K 

where  e^  represents  a  random  error  term  such  that 

E[ek]=0  ;  Var[ek]=(2/N)rk-202; 

Cov[ek,£j]=0   if  k*j    ; 
and  where 

(1)  the   Qk   are   observable    random  variables; 

(2)  the   r^      are  known    constants. 

2 

(3)  9^,  92,  93,  and  a  are  unknown  parameters. 

The  model  can  be  represented  as  a  single  matrix  equation 


Q  =  X9  +  e 


E[e]  =  0 


Cov[ e]=o  V 


(5.0.1) 


where  the  vectors  and  matrices  are 


Qi 
Q332-K 


Q  = 


Q331 
Q332 
Q333 


Q332+K 


x  = 


1 

-K 


0 

-1 

0 

0 

0 

1 

K 
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e  = 


e332-K 


e331 
e332 
e333 


E332+K 


V  =  (1/N) 


1       0 

0    [J_K(16)]-2 


•   •   • 


[J0(16) 


1-2 


[JKU6)]-2 


For  these  specifications  of  the  linear  model,  unbiased  estimators  of  9  and  a 
are  given  by  [8,  page  207] : 


-lv\-lv  "m-^-i 


=    (X'\I~  lX)~  lX'V~  LQ 


»  2 

a  =  J_  Q'[v~1-V"1X(X'V~1X)~1X'V~1]q 
n-3 


(5.0.2) 


Variances  and  covariances  of  the  9's  can  be  estimated  by  substituting  the 
estimator  of  c2  in 


Cov[9]  =  o2(X'V_1X)-1 


(5.0.3) 


The  quantities  listed  above  are  those  needed  for  point  estimation  and 
confidence  intervals  involving  9  and  a2. 


5.1   VOR  Bearing  Angle  Estimator 

According  to  the  reparameterization  of  the  original  nonlinear  model,  the 
relevant  phase  angle  for  VOR  applications  is  9^  .   The  estimator  of  9^  can  be 
obtained  by  algebraic  expansion  of  the  matrix  equation  in  (5.0.2),  and  its 
corresponding  variance  is  the  first  element  of  the  square  matrix  in  (5.0.3). 
The  estimator  of  VOR  bearing  angle  and  its  estimated  variance  are 


where 


K 
h   =  Ql   "  Ck    I   j  Jj2(16)[Q332+j  -  Q332-J 

Vad^]    =   (l/N)(l+CK)a2 


K 


CK  =    [2   I      j2Jj2(16)]-1 

j=l 


(5.1.1) 
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Note  that  the  variance  of  the  bearing  angle  estimator  depends  on  the  selected 
number  of  phase  values  through  C^.   Because  the  Bessel  function  J-j(16) 
approaches  zero  as  j  increases,  it  is  clear  that  C^,  and  hence  Var[  9jJ  , 
approaches  a  lower  bound  as  K  increases.   In  practice,  computational  speed 
and/or  memory  constraints  may  require  that  only  a  few  phase  observations  be 
used  to  estimate  9j  .   Values  of  C -j  listed  in  Table  2  suggest  that  K=10  is 
sufficient  to  minimize  C^  . 


Table  2.   Bessel  Functions  and  Weights 


Jj2(16) 


;2T  2 


JZJ/(16) 


0 

1 

2 

3 

4 

5 

6 

7 

8 

9 

10 

11 

12 

13 

14 

15 

16 

17 

18 

19 

20 

21 

22 

23 

24 


0.0305896861 
0.0081716494 
0.0346699637 
0.0019226029 
0.0410635904 
0.0033031768 
0.0277958044 
0.0333112958 
0.0000492964 
0.0359235032 
0.0425207884 
0.0046542621 
0.0126338128 
0.0560848988 
0.0742215568 
0.0575717228 
0.0314896359 
0.0132170213 
0.0044686657 
0.0012561973 
0.0003002854 
0.0000620785 
0.0000112467 
0.0000018048 
0.0000002588 


0.0000000000 
0.0081716494 
0.1386798547 
0.0173034257 
0.6570174461 
0.0825794204 
1.0006489577 
1.6322534965 
0.0031549718 
2.9098037600 
4.2520788410 
0.5631657107 
1.8192690450 
9.4783478893 
14.5474251288 
12.9536376362 
8.0613467818 
3.8197191653 
1.4478476965 
0.4534872292 
0.1201141759 
0.0273766197 
0.0054433896 
0.0009547344 
0.0001490789 


61.2080 
3.4032 
3.0448 
.5872 
.5360 
.2800 
.1264 
.1264 
.0752 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 
.0240 


However,  the  apparent  gain  from  using  using  only  a  few  phase  variates  is 
balanced  by  a  corresponding  loss  in  precision  for  estimating  a^   .   According 
to  the  specifications  for  equipment  designed  to  provide  the  sample  values  from 
a  VOR  signal,  it  is  likely  that  K=10  will  provide  adequate  precision  for 
estimating  8^  . 

The  procedure  described  above  can  be  used  to  achieve  acceptable  precision 
bounds  using  only  a  few  lines  in  the  phase  spectrum  of  the  observed  VOR 
signal.   An  alternative  method,  which  may  require  fewer  observations,  is  to 
obtain  an  estimate  of  Q\   from  a  few  phase  variables  not  adjacent  to  Q332  •   To 
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illustrate  this  approach,  values  of  j  J.;2 (16)  that  appear  in  CK  are  listed  in 
Table  2.   The  maximum  value  of  this  quantity  occurs  if  j  =  14,  corresponding  to 
phase  observations  at  k=318  and  k=346.   Clearly,  to  minimize  the  variance  of 
an  estimator  of  8^  based  on  nonadjacent  Q^'s,  one  should  add  observations  in 
order  of  decreasing  values  on  j  J  -s  (16).   Thus,  phase  variate  pairs  would  be 
included  in  the  order  j=14, 15, 13, 16,  etc.,  corresponding  to  k=(318,346), 
(317, 347), (319, 345) ,(316,348) ,  etc.   It  is  easy  to  show  that  the  multiplier 
analogous  to  Cj^  is  already  near  the  lower  bound  of  Table  2  after  only  four 
phase  pairs  are  included  to  estimate  9^  .   Hence,  if  estimation  of  a2  is  not 
severely  affected,  a  significant  saving  in  computational  requirements  is 
achieved  using  phase  variables  nonadjacent  to  Q332  •   Assuming  that  N  is 
moderately  large,  the  gain  is  especially  desirable  if  Discrete  Fourier 
Transforms  are  used  to  obtain  the  Q^. 

5.2   Discussion 

It  can  be  argued  that  a  transformation  of  the  nonlinear  model  to  linear 
form  is  unnecessary  because  suitable  nonlinear  least  squares  methods  can  be 
applied  directly  to  the  sample  data.   These  methods  are  iterative  and  require 
initial  estimates  of  the  unknown  phase  angles.   However,  for  VOR  applications, 
software  for  estimating  unknown  angles  will  be  implemented  on  desktop 
computers,  which  will  also  serve  as  controllers  in  VOR  calibration  systems. 
In  this  case,  the  phase  spectrum  transformation  and  subsequent  estimation  of 
phase  angles  using  the  linear  models  approach  is  computationally  efficient  and 
is  to  be  preferred  if  there  are  no  serious  deficiencies  in  the  technique. 

From  a  mathematical  standpoint,  estimation  of  _8  using  the  phase  spectrum 
transformation  depends  on  two  related  assumptions.   First,  it  was  implicitly 
assumed  that  formulas  for  the  mean  and  variance  of  phase  random  variables 
approximate  the  true  mean  and  variance  to  the  extent  that  departures  from  the 
correct  values  are  negligible.   A  second  assumption,  which  requires  further 
study,  concerns  the  adjustment  of  computed  values  of  the  Arctan  function. 
Recall  that  Qk^Qk*"^^71'   where  Q^*  is  an  observable  random  variable  in  the 
interval  (  — tt ,  tt]  .   In  the  derivation  of  the  estimators  of  _9  and  a2,  it  was 
assumed  that  the  integers  (Y^}  can  be  determined  from  the  data.   If  the  Y^'s 
are  not  known,  then  8  and  a2  in  (5.0.2)  are  not  estimators  because  they  are 
not  observable. 

A  computer  simulation  of  the  VOR  signal  with  independent  Gaussian  errors 
was  used  to  determine  if  the  assumptions  described  above  severely  limit 
applicability  of  the  phase  spectrum  transformation.   Results  of  this 
investigation  indicate  that  linearity  of  the  mean  and  the  ability  to  adjust 
the  Arctan  function  depend,  primarily,  on  variability  in  measurement  errors. 
The  technique  was  applied  consistently  for  values  of  0  less  than  .001.   For 
somewhat  greater  values  of  a2,  straightforward  determination  of  the  Y's  is 
usually  successful,  and  indications  are  that  the  method  can  be  refined, 
perhaps  by  developing  a  search  technique  for  the  Y's.   Analyses  were  conducted 
with  N=1024  time  samples. 

It  should  be  emphasized  that  values  of  o^   used  in  computer  simulations 
are  believed  to  be  representative  of  expected  variability  of  measurement 
errors  for  a  VOR  audiofrequency  generator  currently  being  constructed. 
Computations  based  on  simulated  data  were  used  to  affirm  the  mathematical 
results  of  previous  sections  and  are  not  reported  here. 
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To  conclude  this  section  we  remark  that  hardware  specifications  for  VOR 
generators  and  the  method  for  sampling  the  continuous  time  signal  will  togeth- 
er determine  the  accuracy  and  reliability  of  specified  values  for  frequency, 
amplitude,  offset,  and  modulation  index.   Because  properties  of  measurement 
errors,  such  as  stability  and  independence,  can  be  affected  by  hardware  and 
software  specifications,  examination  of  estimated  residuals  for  the  time 
samples  can  be  useful  to  validate  assumptions  about  sampling  errors. 
Estimated  residuals  for  the  VOR  model  are  given  by 

sj  =  Yj  "  yj(i}  »  j=i»---»N    • 

Plots  of  residuals  and/or  tests  for  serial  correlation  can  be  expected  to 
reveal  inconsistent  or  unusual  properties  of  a  particular  VOR  measurement 
system.   Detection  of  a  problem  may  require  redesign  of  the  system  or  a 
modification  of  the  estimation  method  developed  in  this  report. 


6.   ROBUSTNESS  OF  ASSUMPTIONS 


In  this  section,  two  generalizations  of  the  VOR  model  are  examined  to 
understand  the  consequences  if  values  of  some  parameters  assumed  to  be  known 
are  in  error.   To  facilitate  the  discussion,  we  state  a  modification  of  the 
VOR  model  which  is  sufficient  for  the  generalizations  considered  in  this 
section: 


Y(t)   =  v(t;6)  +  s(t;9)  +  e(t) 


where 


and 


v(t;_9)   =  04    cos[2Tr30t  +   Q]_  +   62  ] 

s(t;_6)    =   ct2    cos[2T:9960t  +   16sin[2ir30t  +   62]  ] 

1/9 

Recall  that  for  VOR  applications  we  assumed  that  aj_  =  a2=21'   . 

6. 1   Signal  Offset 

Let  Y(t)  denote  a  signal  with  a^=a2=21/2  .   Suppose  that  instead  of  Y(t) 
we  observe  Y(t)+M  where  M*0.   The  Fourier  transform  of  the  observed  process 
is: 

1/60 
60   /   JY(t)  +  u]exp[i2Tr30kt]  dt 
-1/60 


,  k=0 

,   k*0 

where  S^.  is  the  transform  of  Y(t).   Therefore,  signal  offset  does  not  affect 
the  previous  results  because  S0  was  not  used  tQ  estimate  _9. 
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6.2   Signal  Amplitude 

In  subsection  3.1  we  obtained  the  Fourier  coefficients  of  a  signal  more 
general  than  the  VOR  requires.   The  derivation  of  essential  results  that 
followed  in  no  way  depended  on  the  particular  values  of  the  amplitude 
parameters  a^  and  o^.   If  y(t  ;_9)=v(t  ;_9)+s(t  ;_9)  denotes  a  VOR  signal  with  o^ 
and  oi2  unspecified,  and  s^a^+ib^  denotes  the  corresponding  transform  of 
y(t;9),  then  it  is  easily  shown  that 


ak 


and 


bk  = 


•  aj  Cos[e1  +  e2]  ,  k=i 

-  a2   Jk-332(16)cos[Q3+(k"332)e2]»  k>2 

-o^  sin[61  +  82]  ,  k=l 

-a2   Jk_332(16)sin[63+(k-332)e2],  k>2 


Clearly,  phase  computations  using  tan(q^)=-bjc/a|c,  which  are  fundamental  to  the 
estimation  method,  are  invariant  with  respect  to  particular  values  of  oij  and 
c*2,  even  if  otj^o^  . 
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APPENDIX 


The  appendix  presents  the  derivation  of  Fourier  coefficients  in  equation 
(3.1.1),  and  obtains  the  distribution  of  phase  variables  when  errors  in  the 
VOR  model  are  Gaussian. 


A. 1   Derivation  of  Fourier  Coefficients 

In  this  section  we  determine  the  Fourier  coefficients  for  the  function 

y(t;_6)  =  P  +  v(t;_6)  +  s(t;_9)  (A.l.l) 

where 

v(t;_8)  =  ax  cos^^t  +  9X  +  02  ] 

and 

s(t;_9)  =  a2  cos^irmfx  +  93  +  3  sin[2Trf1  +  92]  ] 

The  following  results  and  trigonometric  identities  are  needed  in  the 
derivation: 


l/2f  f l/2f    ,      j=k*0 

/      cos(2TTfjt)cos(2Trfkt)dt   =<  (A.1.2) 

■l/2f  lO  ,      j^k 


l/2f 

/      cos(2Trf jt)sin(2irfkt)dt   =     0,      all   j,k  (A.1.3) 

•l/2f 


cos((*f3)    =  cos   a  cos    3  -   sin   a  sin   3  (A. 1.4) 

cos   a  cos    3  =   (1/2)  [cos(a-3)   +  cos(oH-3)]  (A.  1.5) 


IT 

/cos(nt   -   z    sin  t)dt   =    TfJn(z)    .  (A.1.6) 

0 


The  Fourier    coefficients    of   y(t;_9)   are   obtained  from  the    complex   integral 

l/2fx 
ak  +  i\  =  2fj_      /     y(t;_e)exp[i2Trf1kt]dt 
-l/2f! 
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where   aj.   and   b^    can   be    represented   by 

ak  =  ak(y)   +  ak^v)    +  ak^s)   1 

\  k=0,±l,±2,'". 

bk  =  bk(y)  +  bk(v)   +  bk(s)    J 

The    components    of   ak  and   b^   represent  Fourier    coefficients    of    respective    terms 
in   (A. 1.1).      Since    p   is   a    constant,    it    is    clear   that 

f    2p      ,        k=0 
ak(M)    = 

I    0         ,        k*0 

and 

bk(y)  0        ,        all  k 

Following  is  a  derivation  of  ak(v)  and  ak(s).   The  corresponding  Fourier  sine 
transforms,  bk(v)  and  ty^Cs),  are  easily  deduced  from  these  results. 

For   non-negative    integer   values   of   k,    we   have 

l/2fx 

ak(v)    =   2f:  /        ax    cos[2Trf1t  +    9^    62  ]cos(2Trf  xkt)dt 

-l/2fj 

l/2f  ! 

=   2f1a1    cos(91   +   02)         /      cos(2irf1t)cos(2TTf1kt)dt 

-l/2fi 

l/2fi 

-    2f1a1sin(e1   +   92)      /        sin(2irf  1t)cos(2uf xkt)dt 
-l/2f1 


r      ax    cos(61   +   92)    ,      k=l 
I     0  ,      k*l 

where  we   have   used   (A. 1.2)    to   (A. 1.4).      Similarly, 


r  -ax    sin(91  +   92)    ,      k=l 
bk(v)    = 

I     0  ,      k*l 

The  Fourier    coefficients   ak(s)   are   given   by 

l/2fx 
ak(s)    =  2fx      /        02    cos[2Trmf1t  +   93  +   3  sin[2irf1t  +    92]  ]cos(2Trf1kt)dt. 
-l/2fx 
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Substituting  u=27rf1t+92,    so   dt=(2irfi)    1t    we   get 

IT 

ak(s)    =   ff-1a2      /   cos[mu   -  m92  +63+3  sin  u]cos(ku  -  k82)   du, 

-Tf 

where   the    limits   -Tf+02<u<Tf+82      can   be   replaced   by   -Tr<u<TT      because   the 
integrand  has    period   2tt.      Furthermore,    using   the   identity   (A.  1.5),    we    obtain 

TT 

ak(s)    =   (2ir)-1a2      /    cos[(m-k)u   -   (m-k)62  +63  +   6   sin  u]du 

-Tf 


+      /    cos[(m+k)u   -   (m+k)62  +63+3   sin  u]du, 

-TT 

and  using    (A. 1.4)   we   get 

Tf 

ak(s)    =   (2Tf)-1a2    cos[03   -   (m-k)82J    J    cos[(m-k)u  +    3  sin  ujdu 

-Tf 

Tf 

-  sin[03   -   (m-k)82J    J    sin[(m-k)u  +   3  sin  ujdu 

-Tf 

Tf 

+   cos [83   -   (m+k)82J    /   cos[(m+k)u  +    3  sin  ujdu 

-Tf 

Tf 

-  sin[83   -   (m+k)62J    J    sin[(m+k)u  +   3  sin  ujdu. 

-Tf 

The  second  and  fourth  integrals  are  zero  because  the  integrands  are  odd 
functions  with  period  2tt.   Since  the  first  and  third  integrands  are  even 
functions,  it  follows  from  (A. 1.6)  that 

ak(s)  =   a2{jm_k(3)cos[83-(m-k)82]  +  Jmhk( 3)cos [e3-(m+k) 82 ] }  . 
Similarly, 

bk(s)  =  -a2{jm_k(3)sin[83-(m-k)82]  -  Jm+k( 3)sin [ 03-(m+k) 82 ] }  . 
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We   have   proved   that    the  Fourier    coefficients   of   y(t;_9)  are 

r  v  ,   k=0 

ak  =  "S     al    cos[Ql   +   e2^    +   a2  ak^s^    »    k=1 
^  a2   ak(s)  ,    k>2 


and 


^  0  ,    k=0 

bk  =^~al    sin[el   +   e2^    ~   a2   bk(s)    >    k=1 


vr«2   bk(s) 
where   a^Cs)   and   b^Cs)   are   defined  above. 


,  k>2 


A. 2  Distribution  of  Phase  Random  Variables 

A  justification  for  using  approximate  formulas  for  the  mean  and  variance 
of  phase  random  variables  is  the  complexity  of  the  exact  distribution  of  the 
Qk* ' s  even  when  errors  are  assumed  to  be  Gaussian.   For  completeness,  the 
distribution  of  Ar ctanf-B^/A^]  is  derived  in  this  section.   The  subscript  is 
dropped  in  the  derivation. 

If  errors  are  Gaussian,  then  A  and  B  are  independent  Gaussian  random 
variables  and 

E[A]=a;  E[B]=b;  a2+b2=r2; 

Var[A]=Var[B]=(2/N)a2. 

Particular  values  of  a  and  b  are  given  by  (4.1.1).   The  joint  distribution  of 
A  and  B  is 

fA)B(u,v)  =  (N/4na2)exp{-(N/4a2)[(u-a)2+(v-b)2]},  -°°<u,v<°°. 

Let   X=A  and  Y=Ar ctan[-B/A] .      Then   because   u=x   and  v=-x   tan  y,    the  Jacobian   of 
the   transformation   is  J=|x|sec  y,   and   the    joint   distribution   of   X  and  Y   is 
given   by 

f(x,y)  =  (N/4a2)  |x|sec2y   exp  {-(N/4a2)  [(x-a)  2+(x   tan  y+b)2]},    -°°<x<«>,-7r/2<y<TT/2. 

Expanding   the   exponent   and    completing  the   square,    we    obtain 

1/9 

f(x,y)    =   (Nsec2y/4iTa2)        exp  {(N/4a2)  [  (a-b   tan  y)2cos2y-r2]  } 


1/9  9 

•  |x|  (Nsec2y/4ira2)        exp  {(-N/4a2)sec2y [x-(a-b   tan  y)cos2y]     } 
=  K(y)«|x|  <f>[x:(a-b   tan  y  )cos2y,  (2/N)  a2cos2y] 
where    <|>[z:5,t    ]    denotes    the   probability  density    function   of   a  Gaussian 
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distribution  with  mean    £  and  variance    t '   .      Integration   of    f(x,y)   over   x  to 
obtain   the   distribution  of  Y   shows    that   f (y )=K(y)E [ |x| ] ,    where   X  is  Gaussian 
with  mean  (a-b  tan  y)cos   y  and  variance  (2/N)o  cos   y. 

If   Z   is   a  Gaussian   random  variable  with  mean    £  and  variance    i ,    then 

0  °° 

E[|Z|]    =  -  /ztf)[z:C,T2]dz  +  Jz<f.[z:?,T2]dz 

-00  0 

00 

=  .-£  +  2/z<f>[z:£,T2]dz. 
0 

Integration   by    parts   gives 

E[|Z|]    =   ?-(l-2$[-C/x:?,x2])   +   (2T2/7T)1/2exp{-?2/2T2}. 

w 
where   $(w:  ?,  t2)=/<J>[z:  £,  T2]dz.      It   follows   that   the  distribution  of 

—  00 

Y=Arctan[-B/A]  is  given  by  f (y)=K(y)E [ |Z | ]  with  £=(a-b  tan  y)cos2y  and 
t  -(2/N)o  cos  y.  Substition  of  these  values  above  gives  the  following 
distribution   for   a   phase   random  variable: 

1   /  9 

f(y)   =  exp(-Nr2/4a2){iT~1+(Ncos2y/4TTa2)        (a-b  tan  y) 

•(l-$[-(Ncos2y/2a2)1/2(a-b  tan  y):?,x2]) 
•exp[(Ncos2y/4a2)(a-b  tan  y)    ]},    -TT/2<y<Tr/2, 
where    £  and    t^  are   defined  above. 
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data  on  the  physical  and  chemical  properties  of  materials,  com- 
piled from  the  world's  literature  and  critically  evaluated. 
Developed  under  a  worldwide  program  coordinated  by  NBS  under 
the  authority  of  the  National  Standard  Data  Act  (Public  Law 
90-396). 


NOTE:  The  principal  publication  outlet  for  the  foregoing  data  is 
the  Journal  of  Physical  and  Chemical  Reference  Data  (JPCRD) 
published  quarterly  for  NBS  by  the  American  Chemical  Society 
(ACS)  and  the  American  Institute  of  Physics  (AIP).  Subscriptions, 
reprints,  and  supplements  available  from  ACS,  1 155  Sixteenth  St., 
NW,  Washington,  DC  20056. 

Building  Science  Series — Disseminates  technical  information 
developed  at  the  Bureau  on  building  materials,  components, 
systems,  and  whole  structures.  The  series  presents  research  results, 
test  methods,  and  performance  criteria  related  to  the  structural  and 
environmental  functions  and  the  durability  and  safety  charac- 
teristics of  building  elements  and  systems. 

Technical  Notes — Studies  or  reports  which  are  complete  in  them- 
selves but  restrictive  in  their  treatment  of  a  subject.  Analogous  to 
monographs  but  not  so  comprehensive  in  scope  or  definitive  in 
treatment  of  the  subject  area.  Often  serve  as  a  vehicle  for  final 
reports  of  work  performed  at  NBS  under  the  sponsorship  of  other 
government  agencies. 

Voluntary  Product  Standards — Developed  under  procedures 
published  by  the  Department  of  Commerce  in  Part  10,  Title  15,  of 
the  Code  of  Federal  Regulations.  The  standards  establish 
nationally  recognized  requirements  for  products,  and  provide  all 
concerned  interests  with  a  basis  for  common  understanding  of  the 
characteristics  of  the  products.  NBS  administers  this  program  as  a 
supplement  to  the  activities  of  the  private  sector  standardizing 
organizations. 

Consumer  Information  Series — Practical  information,  based  on 
NBS  research  and  experience,  covering  areas  of  interest  to  the  con- 
sumer. Easily  understandable  language  and  illustrations  provide 
useful  background  knowledge  for  shopping  in  today's  tech- 
nological marketplace. 

Order  the  above  NBS  publications  from:  Superintendent  of  Docu- 
ments, Government  Printing  Office,   Washington,  DC  20402. 
Order  the  following  NBS  publications — FIPS  and  NBSIR's — from 
the  National  Technical  Information  Services,  Springfield,  VA  22161 . 

Federal    Information    Processing    Standards    Publications    (FIPS 

PUB) — Publications  in  this  series  collectively  constitute  the 
Federal  Information  Processing  Standards  Register.  The  Register 
serves  as  the  official  source  of  information  in  the  Federal  Govern- 
ment regarding  standards  issued  by  NBS  pursuant  to  the  Federal 
Property  and  Administrative  Services  Act  of  1949  as  amended, 
Public  Law  89-306  (79  Stat.  1127),  and  as  implemented  by  Ex- 
ecutive Order  11717(38  FR  12315,  dated  May  11,  1973)  and  Part  6 
of  Title  15  CFR  (Code  of  Federal  Regulations). 

NBS  Interagency  Reports  (NBSIR) — A  special  series  of  interim  or 
final  reports  on  work  performed  by  NBS  for  outside  sponsors 
(both  government  and  non-government).  In  general,  initial  dis- 
tribution is  handled  by  the  sponsor;  public  distribution  is  by  the 
National  Technical  Information  Services,  Springfield,  VA  22161, 
in  paper  copy  or  microfiche  form. 


BIBLIOGRAPHIC  SUBSCRIPTION  SERVICES 


The  following  current-awareness  and  literature-survey  bibliographies 
are  issued  periodically  by  the  Bureau: 

Cryogenic  Data  Center  Current  Awareness  Service.  A  literature  sur- 
vey issued  biweekly.  Annual  subscription:  domestic  $25;  foreign 
$30. 

Liquefied  Natural  Gas.  A  literature  survey  issued  quarterly.  Annual 
subscription:  $20. 


Superconducting  Devices  and  Materials.  A  literature  survey  issued 
quarterly.  Annual  subscription:  $30.  Please  send  subscription  or- 
ders and  remittances  for  the  preceding  bibliographic  services  to  the 
National  Bureau  of  Standards,  Cryogenic  Data  Center  (736) 
Boulder,  CO  80303. 
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